/*
 * Copyright 2020 Leo McCormack
 *
 * Permission to use, copy, modify, and/or distribute this software for any
 * purpose with or without fee is hereby granted, provided that the above
 * copyright notice and this permission notice appear in all copies.
 *
 * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES WITH
 * REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY
 * AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT,
 * INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
 * LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR
 * OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR
 * PERFORMANCE OF THIS SOFTWARE.
 */

/**
 * @file saf_utility_qmf.c
 * @ingroup Utilities
 * @brief An implementation of the complex Quadrature Mirror Filterbank (QMF)
 *        described in [1].
 *
 * @see [1] Herre, J., Purnhagen, H., Breebaart, J., Faller, C., Disch, S.,
 *          Kjörling, K., Schuijers, E., Hilpert, J. and Myburg, F., 2005. The
 *          reference model architecture for MPEG spatial audio coding.
 *
 * @author Leo McCormack
 * @date 14.07.2020
 * @license ISC
 */

#include "saf_utilities.h"
#include "saf_externals.h" 

#define QMF_MAX_HOP_SIZE ( 128 )        /**< Maximum hop size supported */
#define QMF_HYBRID_FILTER_LENGTH ( 13 ) /**< Hybrid-filter length */
#define QMF_NBANDS_2_SUBDIVIDE ( 3 )    /**< Number of QMF bands to subdivide */

/** Prototype filter/window */
static const double __qmf_protofilter[1280] =
{    0.0000000000000, -0.0002762634720, -0.0005525342198, -0.0005571552713, -0.0005617683373, -0.0005282514607, -0.0004947593820, -0.0004911467040, -0.0004875236552, -0.0004884532092, -0.0004893838045, -0.0004967240815, -0.0005040750825, -0.0005133670560, -0.0005226597455, -0.0005346514478, -0.0005466472274, -0.0005572187697, -0.0005677828639, -0.0005774439998, -0.0005870916039, -0.0006001850920, -0.0006132728105, -0.0006222713668, -0.0006312541653, -0.0006426322143, -0.0006540371501, -0.0006658914180, -0.0006777676604, -0.0006859690454, -0.0006941641036, -0.0007049595278, -0.0007157690669, -0.0007206486380, -0.0007255024482, -0.0007348050240, -0.0007440890752, -0.0007465733138, -0.0007490661119, -0.0007586070399, -0.0007681336013, -0.0007703168846, -0.0007724918034, -0.0007779534357, -0.0007834287454, -0.0007807033412, -0.0007779944598, -0.0007791681113, -0.0007803615976, -0.0007802465749, -0.0007801371586, -0.0007779723824, -0.0007758066727, -0.0007694345970, -0.0007630873596, -0.0007580352349, -0.0007530028195, -0.0007424619382, -0.0007319442929, -0.0007267356718, -0.0007215430902, -0.0007066606628, -0.0006917932571, -0.0006784128723, -0.0006650339491, -0.0006495992993, -0.0006341674870, -0.0006143876230, -0.0005946202616, -0.0005755259068, -0.0005564609970, -0.0005355145456, -0.0005145569911, -0.0004875885201, -0.0004606378785, -0.0004350627059, -0.0004095097465, -0.0003798107209, -0.0003501149003, -0.0003199111004, -0.0002897008523, -0.0002497595743, -0.0002098366519, -0.0001772307680, -0.0001446327356, -0.0001031896779, -0.0000617345700, -0.0000241262916, 0.0000134965582, 0.0000614720333, 0.0001094394708, 0.0001568793435, 0.0002042941489, 0.0002496306494, 0.0002949442047, 0.0003487948377, 0.0004026588663, 0.0004567034861, 0.0005107402212, 0.0005673469176, 0.0006239305560, 0.0006848671621, 0.0007458077212, 0.0008033299866, 0.0008608421186, 0.0009247267892, 0.0009886054038, 0.0010567976484, 0.0011250228994, 0.0011914037242, 0.0012577875050, 0.0013240150513, 0.0013902546415, 0.0014672773757, 0.0015443257222, 0.0016155684463, 0.0016868125944, 0.0017608243873, 0.0018348194834, 0.0019094781805, 0.0019841117942, 0.0020651392863, 0.0021461575925, 0.0022239322768, 0.0023017300059, 0.0023821487099, 0.0024625676012, 0.0025413674009, 0.0026201724337, 0.0027036125466, 0.0027870398286, 0.0028669900755, 0.0029469530526, 0.0030297492329, 0.0031125415322, 0.0031932528503, 0.0032739529064, 0.0033579247943, 0.0034418941739, 0.0035213512278, 0.0036008198759, 0.0036806138717, 0.0037603943985, 0.0038405755459, 0.0039207434248, 0.0040013511319, 0.0040819786631, 0.0041542084654, 0.0042264200821, 0.0042997543412, 0.0043730783438, 0.0044470261882, 0.0045209852815, 0.0045908246900, 0.0046606421516, 0.0047269540845, 0.0047932482139, 0.0048535152844, 0.0049137657537, 0.0049765248518, 0.0050393080129, 0.0050900151409, 0.0051407363835, 0.0051934305868, 0.0052461239441, 0.0052966510528, 0.0053471668162, 0.0053834300667, 0.0054196692329, 0.0054536384534, 0.0054875987642, 0.0055175819312, 0.0055475647433, 0.0055706796759, 0.0055937957860, 0.0056079276013, 0.0056220645659, 0.0056337986965, 0.0056455222179, 0.0056422297604, 0.0056389176362, 0.0056327668217, 0.0056266091656, 0.0056091574042, 0.0055917226619, 0.0055660828884, 0.0055404363890, 0.0055079038122, 0.0054753740748, 0.0054296364962, 0.0053838906141, 0.0053277448432, 0.0052715748858, 0.0052048917290, 0.0051382369243, 0.0050611031180, 0.0049839753897, 0.0048974656698, 0.0048109527339, 0.0047074565514, 0.0046039611143, 0.0044920653988, 0.0043801801710, 0.0042526773209, 0.0041251565693, 0.0039854004245, 0.0038456425302, 0.0036928749747, 0.0035401280117, 0.0033746562534, 0.0032091902937, 0.0030269293555, 0.0028446739539, 0.0026477746403, 0.0024508568721, 0.0022391306831, 0.0020274204770, 0.0018029514204, 0.0015784683781, 0.0013343452136, 0.0010902376775, 0.0008367299190, 0.0005832296601, 0.0003054206134, 0.0000275949643, -0.0002594200645, -0.0005464255491, -0.0008516115925, -0.0011568197317, -0.0014803831871, -0.0018039541267, -0.0021433122423, -0.0024826685120, -0.0028380251464, -0.0031933817465, -0.0035667451218, -0.0039401076026, -0.0043311805485, -0.0047222586890, -0.0051279827987, -0.0055337255577, -0.0059564802822, -0.0063792311587, -0.0068204043175, -0.0072615851204, -0.0077207082313, -0.0081798213493, -0.0086561813865, -0.0091325370011, -0.0096237738418, -0.0101150314596, -0.0106232950447, -0.0111315462965, -0.0116582790656, -0.0121849952573, -0.0127284078410, 0.0132718219980, 0.0138311451561, 0.0143904592637, 0.0149655193874, 0.0155405563257, 0.0161365033964, 0.0167324755461, 0.0173379060568, 0.0179433390157, 0.0185652853849, 0.0191872443547, 0.0198202069558, 0.0204531843140, 0.0210999761274, 0.0217467552523, 0.0224073807150, 0.0230680077052, 0.0237420639970, 0.0244160996468, 0.0251018502840, 0.0257875813769, 0.0264867637542, 0.0271859403015, 0.0278965786129, 0.0286072223055, 0.0293287412129, 0.0300502657578, 0.0307760035778, 0.0315017677644, 0.0322385911186, 0.0329754179605, 0.0337187479017, 0.0344620871585, 0.0352159270599, 0.0359697555143, 0.0367255220455, 0.0374812763540, 0.0382433168854, 0.0390053697576, 0.0397701421370, 0.0405349210370, 0.0412999167409, 0.0420649194101, 0.0428373278803, 0.0436097603402, 0.0443793037143, 0.0451488448202, 0.0459165741975, 0.0466843062700, 0.0474504368762, 0.0482165646417, 0.0489775642920, 0.0497385823153, 0.0504970964492, 0.0512556073040, 0.0520093457341, 0.0527630662304, 0.0535041667383, 0.0542452690526, 0.0549813220428, 0.0557173653993, 0.0564395132831, 0.0571616382627, 0.0578766011820, 0.0585915643933, 0.0592876642498, 0.0599837452212, 0.0606646387769, 0.0613455131702, 0.0620156474286, 0.0626857856455, 0.0633286936041, 0.0639715979264, 0.0645981487145, 0.0652247009326, 0.0658307267839, 0.0664367611136, 0.0670221717462, 0.0676076079753, 0.0681559891452, 0.0687043838422, 0.0692337073129, 0.0697630166725, 0.0702629549865, 0.0707628759478, 0.0712315634186, 0.0717002581741, 0.0721342576042, 0.0725682604314, 0.0729651377816, 0.0733620227657, 0.0737312022298, 0.0741003643894, 0.0744228140298, 0.0747452534776, 0.0750294863750, 0.0753137309752, 0.0755572780429, 0.0758008300938, 0.0760000434804, 0.0761992513763, 0.0763492295011, 0.0764992119303, 0.0766042744439, 0.0767093398745, 0.0767633751145, 0.0768173954582, 0.0768201997673, 0.0768229914942, 0.0767717417580, 0.0767204942750, 0.0766127798195, 0.0765050741722, 0.0763399517187, 0.0761748271048, 0.0759527012645, 0.0757305674703, 0.0754440943880, 0.0751576275372, 0.0748120389752, 0.0744664374116, 0.0740535275534, 0.0736406055728, 0.0731590374596, 0.0726774627175, 0.0721300597492, 0.0715826345203, 0.0709679700483, 0.0703533154025, 0.0696598613057, 0.0689663921001, 0.0682094555550, 0.0674525095060, 0.0666107857601, 0.0657690713457, 0.0648567698624, 0.0639444739640, 0.0629523756575, 0.0619602868848, 0.0608884629011, 0.0598166659796, 0.0586659638617, 0.0575152710554, 0.0562806383895, 0.0550459966795, 0.0537276876358, 0.0524093819433, 0.0510036276363, 0.0495978728266, 0.0481141015474, 0.0466303207525, 0.0450536063583, 0.0434768684654, 0.0418113470397, 0.0401458200425, 0.0383938192996, 0.0366418135282, 0.0348001051306, 0.0329583954927, 0.0310203937989, 0.0290824047877, 0.0270565758359, 0.0250307595445, 0.0229152370232, 0.0207997071989, 0.0185849092217, 0.0163701328561, 0.0140662590890, 0.0117623737968, 0.0093630256607, 0.0069636910976, 0.0044701237103, 0.0019765626627, -0.0006160718469, -0.0032086847081, -0.0058899293948, -0.0085711743427, -0.0113500234728, -0.0141288758122, -0.0170061403777, -0.0198834128941, -0.0228530634389, -0.0258227375650, -0.0288879214102, -0.0319531234562, -0.0351153841344, -0.0382776482545, -0.0415291715265, -0.0447806909816, -0.0481305485026, -0.0514804188177, -0.0549254768665, -0.0583705333511, -0.0619057688888, -0.0654409770886, -0.0690676667901, -0.0726943349420, -0.0764158049026, -0.0801372926756, -0.0839460228550, -0.0877547532403, -0.0916540369715, -0.0955533319905, -0.0995431447659, -0.1035329487993, -0.1076078309406, -0.1116826856452, -0.1158452448469, -0.1200078053651, -0.1242540480203, -0.1285002933462, -0.1328277312573, -0.1371551768885, -0.1415659077401, -0.1459766566189, -0.1504686708336, -0.1549607093374, -0.1595283025664, -0.1640958787703, -0.1687383503701, -0.1733808142610, -0.1780990400378, -0.1828172523734, -0.1876069690506, -0.1923966770678, -0.1972608479266, -0.2021250183703, -0.2070492992006, -0.2119735831914, -0.2169694342266, -0.2219652633384, -0.2270171791500, -0.2320690812132, -0.2371853956606, -0.2423016905529, -0.2474748586816, -0.2526480256489, -0.2578766747326, -0.2631053356501, -0.2683843618811, -0.2736634049849, -0.2789924190238, -0.2843214254088, -0.2896965495685, -0.2950716774625, -0.3004907673195, -0.3059098614390, -0.3113688729716, -0.3168278894291, -0.3223196220800, -0.3278113723324, -0.3333418217172, -0.3388722635717, -0.3444318350435, -0.3499914122649, -0.3555752006991, 0.3611589994542, 0.3667692644921, 0.3723795603770, 0.3780072824452, 0.3836350041146, 0.3892780834412, 0.3949211684032, 0.4005764782039, 0.4062317776389, 0.4119007215285, 0.4175696873908, 0.4232408506554, 0.4289119900372, 0.4345836835439, 0.4402553743445, 0.4459275135166, 0.4515996558217, 0.4572652216857, 0.4629308042606, 0.4685880677215, 0.4742453243888, 0.4798853201396, 0.4855253145854, 0.4911480648994, 0.4967708227669, 0.5023762846177, 0.5079817400724, 0.5135526142353, 0.5191235000340, 0.5246737972393, 0.5302240884908, 0.5357397180600, 0.5412553475411, 0.5467302419087, 0.5522051263104, 0.5576420296770, 0.5630789053140, 0.5684656541449, 0.5738524091594, 0.5791963774384, 0.5845403276635, 0.5898263122286, 0.5951123147636, 0.6003453352359, 0.6055783444773, 0.6107446766601, 0.6159109984769, 0.6210176335538, 0.6261242644607, 0.6311611315428, 0.6361980097723, 0.6411624903953, 0.6461269646247, 0.6510142937142, 0.6559016399581, 0.6607078128371, 0.6655139870093, 0.6702401473311, 0.6749663289682, 0.6796008198535, 0.6842353261119, 0.6887817750021, 0.6933282316192, 0.6977835605066, 0.7022388743098, 0.7065899657830, 0.7109410429680, 0.7151936458078, 0.7194462550033, 0.7235955759236, 0.7277448954474, 0.7317830379552, 0.7358211812312, 0.7397519890246, 0.7436827956087, 0.7474982580821, 0.7513137453267, 0.7550109081443, 0.7587080846807, 0.7622877890950, 0.7658674845073, 0.7693227784089, 0.7727780946395, 0.7761034130896, 0.7794287461480, 0.7826320306499, 0.7858353214217, 0.7889044522646, 0.7919735858404, 0.7949101196306, 0.7978466399825, 0.8006476091170, 0.8034485723241, 0.8061090405782, 0.8087694965052, 0.8112943189503, 0.8138191301784, 0.8161983628005, 0.8185776070011, 0.8208098026271, 0.8230419833267, 0.8251347532814, 0.8272275373755, 0.8291656850689, 0.8311038496718, 0.8328987918843, 0.8346937360492, 0.8363327361895, 0.8379717266217, 0.8394629333539, 0.8409541296434, 0.8422889785528, 0.8436238352155, 0.8448028351663, 0.8459818486135, 0.8470067050139, 0.8480315835964, 0.8489060396447, 0.8497805277784, 0.8504888339266, 0.8511971453674, 0.8517509372895, 0.8523047120534, 0.8527033945043, 0.8531021044567, 0.8533370702489, 0.8535720597213, 0.8536553025990, 0.8537385506405, 0.8536553025990, 0.8535720597213, 0.8533370702489, 0.8531021044567, 0.8527033945043, 0.8523047120534, 0.8517509372895, 0.8511971453674, 0.8504888339266, 0.8497805277784, 0.8489060396447, 0.8480315835964, 0.8470067050139, 0.8459818486135, 0.8448028351663, 0.8436238352155, 0.8422889785528, 0.8409541296434, 0.8394629333539, 0.8379717266217, 0.8363327361895, 0.8346937360492, 0.8328987918843, 0.8311038496718, 0.8291656850689, 0.8272275373755, 0.8251347532814, 0.8230419833267, 0.8208098026271, 0.8185776070011, 0.8161983628005, 0.8138191301784, 0.8112943189503, 0.8087694965052, 0.8061090405782, 0.8034485723241, 0.8006476091170, 0.7978466399825, 0.7949101196306, 0.7919735858404, 0.7889044522646, 0.7858353214217, 0.7826320306499, 0.7794287461480, 0.7761034130896, 0.7727780946395, 0.7693227784089, 0.7658674845073, 0.7622877890950, 0.7587080846807, 0.7550109081443, 0.7513137453267, 0.7474982580821, 0.7436827956087, 0.7397519890246, 0.7358211812312, 0.7317830379552, 0.7277448954474, 0.7235955759236, 0.7194462550033, 0.7151936458078, 0.7109410429680, 0.7065899657830, 0.7022388743098, 0.6977835605066, 0.6933282316192, 0.6887817750021, 0.6842353261119, 0.6796008198535, 0.6749663289682, 0.6702401473311, 0.6655139870093, 0.6607078128371, 0.6559016399581, 0.6510142937142, 0.6461269646247, 0.6411624903953, 0.6361980097723, 0.6311611315428, 0.6261242644607, 0.6210176335538, 0.6159109984769, 0.6107446766601, 0.6055783444773, 0.6003453352359, 0.5951123147636, 0.5898263122286, 0.5845403276635, 0.5791963774384, 0.5738524091594, 0.5684656541449, 0.5630789053140, 0.5576420296770, 0.5522051263104, 0.5467302419087, 0.5412553475411, 0.5357397180600, 0.5302240884908, 0.5246737972393, 0.5191235000340, 0.5135526142353, 0.5079817400724, 0.5023762846177, 0.4967708227669, 0.4911480648994, 0.4855253145854, 0.4798853201396, 0.4742453243888, 0.4685880677215, 0.4629308042606, 0.4572652216857, 0.4515996558217, 0.4459275135166, 0.4402553743445, 0.4345836835439, 0.4289119900372, 0.4232408506554, 0.4175696873908, 0.4119007215285, 0.4062317776389, 0.4005764782039, 0.3949211684032, 0.3892780834412, 0.3836350041146, 0.3780072824452, 0.3723795603770, 0.3667692644921, -0.3611589811729, -0.3555752006991, -0.3499914122649, -0.3444318350435, -0.3388722635717, -0.3333418217172, -0.3278113723324, -0.3223196220800, -0.3168278894291, -0.3113688729716, -0.3059098614390, -0.3004907673195, -0.2950716774625, -0.2896965495685, -0.2843214254088, -0.2789924190238, -0.2736634049849, -0.2683843618811, -0.2631053356501, -0.2578766747326, -0.2526480256489, -0.2474748586816, -0.2423016905529, -0.2371853956606, -0.2320690812132, -0.2270171791500, -0.2219652633384, -0.2169694342266, -0.2119735831914, -0.2070492992006, -0.2021250183703, -0.1972608479266, -0.1923966770678, -0.1876069690506, -0.1828172523734, -0.1780990400378, -0.1733808142610, -0.1687383503701, -0.1640958787703, -0.1595283025664, -0.1549607093374, -0.1504686708336, -0.1459766566189, -0.1415659077401, -0.1371551768885, -0.1328277312573, -0.1285002933462, -0.1242540480203, -0.1200078053651, -0.1158452448469, -0.1116826856452, -0.1076078309406, -0.1035329487993, -0.0995431447659, -0.0955533319905, -0.0916540369715, -0.0877547532403, -0.0839460228550, -0.0801372926756, -0.0764158049026, -0.0726943349420, -0.0690676667901, -0.0654409770886, -0.0619057688888, -0.0583705333511, -0.0549254768665, -0.0514804188177, -0.0481305485026, -0.0447806909816, -0.0415291715265, -0.0382776482545, -0.0351153841344, -0.0319531234562, -0.0288879214102, -0.0258227375650, -0.0228530634389, -0.0198834128941, -0.0170061403777, -0.0141288758122, -0.0113500234728, -0.0085711743427, -0.0058899293948, -0.0032086847081, -0.0006160718469, 0.0019765626627, 0.0044701237103, 0.0069636910976, 0.0093630256607, 0.0117623737968, 0.0140662590890, 0.0163701328561, 0.0185849092217, 0.0207997071989, 0.0229152370232, 0.0250307595445, 0.0270565758359, 0.0290824047877, 0.0310203937989, 0.0329583954927, 0.0348001051306, 0.0366418135282, 0.0383938192996, 0.0401458200425, 0.0418113470397, 0.0434768684654, 0.0450536063583, 0.0466303207525, 0.0481141015474, 0.0495978728266, 0.0510036276363, 0.0524093819433, 0.0537276876358, 0.0550459966795, 0.0562806383895, 0.0575152710554, 0.0586659638617, 0.0598166659796, 0.0608884629011, 0.0619602868848, 0.0629523756575, 0.0639444739640, 0.0648567698624, 0.0657690713457, 0.0666107857601, 0.0674525095060, 0.0682094555550, 0.0689663921001, 0.0696598613057, 0.0703533154025, 0.0709679700483, 0.0715826345203, 0.0721300597492, 0.0726774627175, 0.0731590374596, 0.0736406055728, 0.0740535275534, 0.0744664374116, 0.0748120389752, 0.0751576275372, 0.0754440943880, 0.0757305674703, 0.0759527012645, 0.0761748271048, 0.0763399517187, 0.0765050741722, 0.0766127798195, 0.0767204942750, 0.0767717417580, 0.0768229914942, 0.0768201997673, 0.0768173954582, 0.0767633751145, 0.0767093398745, 0.0766042744439, 0.0764992119303, 0.0763492295011, 0.0761992513763, 0.0760000434804, 0.0758008300938, 0.0755572780429, 0.0753137309752, 0.0750294863750, 0.0747452534776, 0.0744228140298, 0.0741003643894, 0.0737312022298, 0.0733620227657, 0.0729651377816, 0.0725682604314, 0.0721342576042, 0.0717002581741, 0.0712315634186, 0.0707628759478, 0.0702629549865, 0.0697630166725, 0.0692337073129, 0.0687043838422, 0.0681559891452, 0.0676076079753, 0.0670221717462, 0.0664367611136, 0.0658307267839, 0.0652247009326, 0.0645981487145, 0.0639715979264, 0.0633286936041, 0.0626857856455, 0.0620156474286, 0.0613455131702, 0.0606646387769, 0.0599837452212, 0.0592876642498, 0.0585915643933, 0.0578766011820, 0.0571616382627, 0.0564395132831, 0.0557173653993, 0.0549813220428, 0.0542452690526, 0.0535041667383, 0.0527630662304, 0.0520093457341, 0.0512556073040, 0.0504970964492, 0.0497385823153, 0.0489775642920, 0.0482165646417, 0.0474504368762, 0.0466843062700, 0.0459165741975, 0.0451488448202, 0.0443793037143, 0.0436097603402, 0.0428373278803, 0.0420649194101, 0.0412999167409, 0.0405349210370, 0.0397701421370, 0.0390053697576, 0.0382433168854, 0.0374812763540, 0.0367255220455, 0.0359697555143, 0.0352159270599, 0.0344620871585, 0.0337187479017, 0.0329754179605, 0.0322385911186, 0.0315017677644, 0.0307760035778, 0.0300502657578, 0.0293287412129, 0.0286072223055, 0.0278965786129, 0.0271859403015, 0.0264867637542, 0.0257875813769, 0.0251018502840, 0.0244160996468, 0.0237420639970, 0.0230680077052, 0.0224073807150, 0.0217467552523, 0.0210999761274, 0.0204531843140, 0.0198202069558, 0.0191872443547, 0.0185652853849, 0.0179433390157, 0.0173379060568, 0.0167324755461, 0.0161365033964, 0.0155405563257, 0.0149655193874, 0.0143904592637, 0.0138311451561, -0.0132718220090, -0.0127284078410, -0.0121849952573, -0.0116582790656, -0.0111315462965, -0.0106232950447, -0.0101150314596, -0.0096237738418, -0.0091325370011, -0.0086561813865, -0.0081798213493, -0.0077207082313, -0.0072615851204, -0.0068204043175, -0.0063792311587, -0.0059564802822, -0.0055337255577, -0.0051279827987, -0.0047222586890, -0.0043311805485, -0.0039401076026, -0.0035667451218, -0.0031933817465, -0.0028380251464, -0.0024826685120, -0.0021433122423, -0.0018039541267, -0.0014803831871, -0.0011568197317, -0.0008516115925, -0.0005464255491, -0.0002594200645, 0.0000275949643, 0.0003054206134, 0.0005832296601, 0.0008367299190, 0.0010902376775, 0.0013343452136, 0.0015784683781, 0.0018029514204, 0.0020274204770, 0.0022391306831, 0.0024508568721, 0.0026477746403, 0.0028446739539, 0.0030269293555, 0.0032091902937, 0.0033746562534, 0.0035401280117, 0.0036928749747, 0.0038456425302, 0.0039854004245, 0.0041251565693, 0.0042526773209, 0.0043801801710, 0.0044920653988, 0.0046039611143, 0.0047074565514, 0.0048109527339, 0.0048974656698, 0.0049839753897, 0.0050611031180, 0.0051382369243, 0.0052048917290, 0.0052715748858, 0.0053277448432, 0.0053838906141, 0.0054296364962, 0.0054753740748, 0.0055079038122, 0.0055404363890, 0.0055660828884, 0.0055917226619, 0.0056091574042, 0.0056266091656, 0.0056327668217, 0.0056389176362, 0.0056422297604, 0.0056455222179, 0.0056337986965, 0.0056220645659, 0.0056079276013, 0.0055937957860, 0.0055706796759, 0.0055475647433, 0.0055175819312, 0.0054875987642, 0.0054536384534, 0.0054196692329, 0.0053834300667, 0.0053471668162, 0.0052966510528, 0.0052461239441, 0.0051934305868, 0.0051407363835, 0.0050900151409, 0.0050393080129, 0.0049765248518, 0.0049137657537, 0.0048535152844, 0.0047932482139, 0.0047269540845, 0.0046606421516, 0.0045908246900, 0.0045209852815, 0.0044470261882, 0.0043730783438, 0.0042997543412, 0.0042264200821, 0.0041542084654, 0.0040819786631, 0.0040013511319, 0.0039207434248, 0.0038405755459, 0.0037603943985, 0.0036806138717, 0.0036008198759, 0.0035213512278, 0.0034418941739, 0.0033579247943, 0.0032739529064, 0.0031932528503, 0.0031125415322, 0.0030297492329, 0.0029469530526, 0.0028669900755, 0.0027870398286, 0.0027036125466, 0.0026201724337, 0.0025413674009, 0.0024625676012, 0.0023821487099, 0.0023017300059, 0.0022239322768, 0.0021461575925, 0.0020651392863, 0.0019841117942, 0.0019094781805, 0.0018348194834, 0.0017608243873, 0.0016868125944, 0.0016155684463, 0.0015443257222, 0.0014672773757, 0.0013902546415, 0.0013240150513, 0.0012577875050, 0.0011914037242, 0.0011250228994, 0.0010567976484, 0.0009886054038, 0.0009247267892, 0.0008608421186, 0.0008033299866, 0.0007458077212, 0.0006848671621, 0.0006239305560, 0.0005673469176, 0.0005107402212, 0.0004567034861, 0.0004026588663, 0.0003487948377, 0.0002949442047, 0.0002496306494, 0.0002042941489, 0.0001568793435, 0.0001094394708, 0.0000614720333, 0.0000134965582, -0.0000241262916, -0.0000617345700, -0.0001031896779, -0.0001446327356, -0.0001772307680, -0.0002098366519, -0.0002497595743, -0.0002897008523, -0.0003199111004, -0.0003501149003, -0.0003798107209, -0.0004095097465, -0.0004350627059, -0.0004606378785, -0.0004875885201, -0.0005145569911, -0.0005355145456, -0.0005564609970, -0.0005755259068, -0.0005946202616, -0.0006143876230, -0.0006341674870, -0.0006495992993, -0.0006650339491, -0.0006784128723, -0.0006917932571, -0.0007066606628, -0.0007215430902, -0.0007267356718, -0.0007319442929, -0.0007424619382, -0.0007530028195, -0.0007580352349, -0.0007630873596, -0.0007694345970, -0.0007758066727, -0.0007779723824, -0.0007801371586, -0.0007802465749, -0.0007803615976, -0.0007791681113, -0.0007779944598, -0.0007807033412, -0.0007834287454, -0.0007779534357, -0.0007724918034, -0.0007703168846, -0.0007681336013, -0.0007586070399, -0.0007490661119, -0.0007465733138, -0.0007440890752, -0.0007348050240, -0.0007255024482, -0.0007206486380, -0.0007157690669, -0.0007049595278, -0.0006941641036, -0.0006859690454, -0.0006777676604, -0.0006658914180, -0.0006540371501, -0.0006426322143, -0.0006312541653, -0.0006222713668, -0.0006132728105, -0.0006001850920, -0.0005870916039, -0.0005774439998, -0.0005677828639, -0.0005572187697, -0.0005466472274, -0.0005346514478, -0.0005226597455, -0.0005133670560, -0.0005040750825, -0.0004967240815, -0.0004893838045, -0.0004884532092, -0.0004875236552, -0.0004911467040, -0.0004947593820, -0.0005282514607, -0.0005617683373, -0.0005571552713, -0.0005525342198, -0.0002762634720    };

/** Coefficients to subdivide the lowest QMF band into 8 subbands:
 *  ~  0.125*kaiser(13, 4.6) */
static const double __qmf_fb8bandCoeffs[13] =
{    0.0074608293425, 0.0227042095837, 0.0454686594059, 0.0726611394476, 0.0988510858112, 0.1179371056821, 0.1250000000000, 0.1179371056821, 0.0988510858112, 0.0726611394476, 0.0454686594059, 0.0227042095837, 0.0074608293425    };

/** Coefficients to subdivide the 2nd & 3rd lowest QMF bands into 4 subbands */
static const double __qmf_fb4bandCoeffs[13] =
{    0.0000000000000, 0.0189948752516, 0.0000000000000, -0.0729313916667, 0.0000000000000, 0.3059663054427, 0.5000000000000, 0.3059663054427, 0.0000000000000, -0.0729313916667, 0.0000000000000, 0.0189948752516, 0.0000000000000    };

/**
 * Matrix for converting the centre frequencies of the first 3 QMF bands into
 * the centre frequencies for the 10 hybrid bands. */
static const float __qmf2hybCentreFreq[10][QMF_NBANDS_2_SUBDIVIDE] =
{ { 0.1013f, 0.0f,    0.0f},
  { 0.2027f, 0.0f,    0.0f},
  { 0.4054f, 0.0f,    0.0f},
  { 0.8108f, 0.0f,    0.0f},
  { 1.2533f, 0.0f,    0.0f},
  { 1.7227f, 0.0f,    0.0f},
  { 0.0f,    0.9039f, 0.0f},
  { 0.0f,    1.1228f, 0.0f},
  { 0.0f,    0.0f,    0.9424f},
  { 0.0f,    0.0f,    1.0672f} };

/** Passes input time-domain data through the QMF filterbank. */
static void qmfAnalyse
(
    float* inTD/* nSamplesTD x nCH */,
    int nSamplesTD,
    int nCH,
    int hopSize,
    int hybridmode,
    float_complex* outTF /* out_nBands x nTimeslots x nCH */
)
{
    int t, ch, sample, band;
    void* hQMF;
    float_complex*** FrameTF;
    float** tempFrameTD;
    int nTimeSlots, nBands;

    nBands = hopSize + (hybridmode ? 7 : 0);
    nTimeSlots = nSamplesTD/hopSize;

    /* allocate memory */
    qmf_create(&(hQMF), nCH, 1, hopSize, hybridmode, QMF_TIME_CH_BANDS);
    FrameTF = (float_complex***)malloc3d(nTimeSlots, nCH, nBands, sizeof(float_complex));
    tempFrameTD = (float**)calloc2d(nCH, nTimeSlots*hopSize, sizeof(float));

    /* perform TF transform */
    for(ch=0; ch<nCH; ch++)
        for(sample=0; sample<nSamplesTD; sample++)
            tempFrameTD[ch][sample] = inTD[sample* nCH + ch];
    qmf_analysis(hQMF, tempFrameTD, nTimeSlots*hopSize, FrameTF);

    /* save result to output */
    for (band = 0; band < nBands; band++)
        for (t = 0; t < nTimeSlots; t++)
            for (ch = 0; ch < nCH; ch++)
                outTF[band * nTimeSlots * nCH + t * nCH + ch] = FrameTF[t][ch][band];

    /* clean-up */
    qmf_destroy(&hQMF);
    free(FrameTF);
    free(tempFrameTD);
}

/** Data structure for the complex-QMF filterbank */
typedef struct _qmf_data {
    int hopsize;                      /**< Hop size in samples */
    int hybridmode;                   /**< 1: hybrid filtering enabled; 0: disabled */
    int nCHin;                        /**< Number of input channels */
    int nCHout;                       /**< Number of output channels*/
    int nBands;                       /**< Number of frequency bands */
    int procDelay;                    /**< Processing delay in samples */
    QMF_FDDATA_FORMAT format;         /**< see #QMF_FDDATA_FORMAT */ 

    /* QMF Analysis and Synthesis filters */
    float_complex** h_a;
    float** h_s_real;
    float** h_s_imag;

    /* Prototype window */
    float* h_p;

    /* For run-time */
    float** buffer_ana;
    float** buffer_syn;
    float* buffer_win;
    float* win_sum;
    float* win_sum_cmplx_dummy; /* treated as complex data type (interleaved with zeros for imag parts) */
    float_complex* qmfTF_frame;
    float* qmfTF_frame_tmp; /* used when taking the real/imag parts */
    float* tmp_real_frame;
    float* tmp_imag_frame;

    /* For hybrid filtering */
    float_complex fb8bandCoeffs[8][QMF_HYBRID_FILTER_LENGTH];
    float_complex fb4bandCoeffs[2][QMF_HYBRID_FILTER_LENGTH];
    float_complex*** hybBuffer;
    float_complex*** qmfDelayBuffer;
    float_complex* hybQmfTF_frame;

}qmf_data;


void qmf_create
(
    void ** const phQMF,
    int nCHin,
    int nCHout,
    int hopsize,
    int hybridmode,
    QMF_FDDATA_FORMAT format
)
{
    *phQMF = malloc1d(sizeof(qmf_data));
    qmf_data *h = (qmf_data*)(*phQMF);
    int i,j,K,N,dsFactor;
    float scale, eq;
    float* k_tmp, *n_tmp;

    saf_assert(hopsize==4 || hopsize==8 || hopsize==16 || hopsize==32 || hopsize==64 || hopsize==128, "Unsupported hopsize");

    h->nCHin = nCHin;
    h->nCHout = nCHout;
    h->hopsize = hopsize;
    h->hybridmode = hybridmode;
    h->nBands = hybridmode ? hopsize+7 : hopsize; /* hybrid mode incurs an additional 7 bands */
    h->format = format;
    
    K = hopsize;
    N = 2*hopsize;
    k_tmp = malloc1d(K*sizeof(float));
    n_tmp = malloc1d(N*sizeof(float));

    /* QMF Analysis filters */
    h->h_a = (float_complex**)malloc2d(K, N, sizeof(float_complex));
    scale = (float)QMF_MAX_HOP_SIZE / (2.0f*(float)hopsize); /* (Used to balance the levels between different hopsizes) */
    for(i=0; i<K; i++)
        k_tmp[i] = SAF_PI/2.0f/(float)K * ((float)i+0.5f);
    for(i=0; i<N; i++)
        n_tmp[i] = 2.0f*(float)i - 2.0f*(float)K/(float)QMF_MAX_HOP_SIZE;
    for(i=0; i<K; i++)
        for(j=0; j<N; j++)
            h->h_a[i][j] = crmulf(cexpf(cmplxf(0.0f, k_tmp[i]*n_tmp[j])), scale);

    /* QMF Synthesis filters */
    h->h_s_real = (float**)malloc2d(N, K, sizeof(float));
    h->h_s_imag = (float**)malloc2d(N, K, sizeof(float));
    scale = 2.0f / QMF_MAX_HOP_SIZE; /* (Used to balance the levels between different hopsizes) */
    for(i=0; i<N; i++)
        n_tmp[i] = 2.0f*(float)i - (2.0f*(float)QMF_MAX_HOP_SIZE-1.0f)*(float)K/((float)QMF_MAX_HOP_SIZE/2.0f);
    for(i=0; i<N; i++){
        for(j=0; j<K; j++){
            h->h_s_real[i][j] = scale * cosf(k_tmp[j]*n_tmp[i]);
            h->h_s_imag[i][j] = scale * sinf(k_tmp[j]*n_tmp[i]);
        }
    }

    /* Prototype filter */
    h->h_p = malloc1d(10*hopsize*sizeof(float));
    if(hopsize<=QMF_MAX_HOP_SIZE)
        for(i=0,j=0; i<10*QMF_MAX_HOP_SIZE; i+=QMF_MAX_HOP_SIZE/hopsize, j++)
            h->h_p[j] = (float)__qmf_protofilter[i];
    else{ /* Borrow the one from afSTFT: */
        eq = 2.0f/sqrtf(5.487604141f);
        dsFactor = 1024/hopsize;
        for(i=0; i<10*hopsize; i++)
            h->h_p[i] = __afSTFT_protoFilter1024[i*dsFactor]*eq;
    }

    /* Run-time buffers */
    h->buffer_ana = (float**)malloc1d(nCHin*sizeof(float*)); /* actually faster to go with non-contiguous allocated memory here due to memmove... */
    for(i=0; i<nCHin; i++)
        h->buffer_ana[i] = calloc1d(hopsize * 10, sizeof(float));
    h->buffer_syn = (float**)malloc1d(nCHout*sizeof(float*));
    for(i=0; i<nCHout; i++)
        h->buffer_syn[i] = calloc1d(hopsize * 20, sizeof(float));
    h->buffer_win = malloc1d(hopsize * 10 * sizeof(float));
    h->win_sum = malloc1d(hopsize * 2 * sizeof(float));
    h->win_sum_cmplx_dummy = calloc1d(hopsize * 4, sizeof(float)); /* ca */
    h->qmfTF_frame = malloc1d(hopsize * sizeof(float_complex));
    h->qmfTF_frame_tmp = malloc1d(hopsize * sizeof(float));
    h->tmp_real_frame = malloc1d(hopsize * 2 * sizeof(float));
    h->tmp_imag_frame = malloc1d(hopsize * 2 * sizeof(float));

    /* Init hybrid filtering coefficients: */
    if(hybridmode){
        /* Coefficients to subdivide the lowest QMF band into 8 subbands. The
         * 1st, 2nd, 7th and 8th subbands then become the 1st, 2nd, 5th and 6th
         * hybrid bands directly. Hybrid bands 3 and 4 are then obtained by
         * summing subbands 3 and 6, and 4 and 5, respectively. Therefore, this
         * raises the total number of bands by 5.
         * (Note: __fb8bandCoeffs = 0.125*kaiser(13, 4.6) ) */
        for(i=0; i<8; i++)
            for(j=0; j<QMF_HYBRID_FILTER_LENGTH; j++)
                h->fb8bandCoeffs[i][j] = crmulf(cexpf(cmplxf(0.0f, -1.0f*SAF_PI*((float)j-((float)QMF_HYBRID_FILTER_LENGTH-1.0f)/2.0f)/8.0f*(1.0f+2.0f*(float)i))),
                                                (float)__qmf_fb8bandCoeffs[j]);

        /* Coefficients to subdivide the 2nd and 3rd lowest QMF bands into 4
         * subbands each. The 1st and 4th subdivided bands, and 2nd and 3rd
         * subdivided bands are then summed to obtain the 2 hybrid bands for
         * each QMF band . Therefore, this raises the total number of bands by
         * 2 */
        for(i=0; i<2; i++)
            for(j=0; j<QMF_HYBRID_FILTER_LENGTH; j++)
                h->fb4bandCoeffs[i][j] = cmplxf((float)__qmf_fb4bandCoeffs[j] *
                                                cosf(2.0f*SAF_PI*(float)i*((float)j-(((float)QMF_HYBRID_FILTER_LENGTH-1.0f)/2.0f))/2.0f), 0.0f);

        /* For run-time */
        h->qmfDelayBuffer = (float_complex***)calloc3d(nCHin, hopsize-QMF_NBANDS_2_SUBDIVIDE, (QMF_HYBRID_FILTER_LENGTH-1)/2 + 1, sizeof(float_complex)); /* ca */
        h->hybBuffer = (float_complex***)calloc3d(nCHin, QMF_NBANDS_2_SUBDIVIDE, QMF_HYBRID_FILTER_LENGTH, sizeof(float_complex));
        h->hybQmfTF_frame = malloc1d(h->nBands * sizeof(float_complex));

        /* Processing delay */
        h->procDelay = hopsize*15+1;
    }
    else{
        /* Processing delay */
        h->procDelay = hopsize*9+1;
    }

    /* clean-up */
    free(k_tmp);
    free(n_tmp);
}

void qmf_destroy
(
    void ** const phQMF
)
{
    qmf_data *h = (qmf_data*)(*phQMF);
    int i;

    if(h!=NULL){
        /* QMF Analysis and Synthesis filters */
        free(h->h_a);
        free(h->h_s_real);
        free(h->h_s_imag);

        /* Prototype window */
        free(h->h_p);

        /* For run-time */
        for(i=0; i<h->nCHin; i++)
            free(h->buffer_ana[i]);
        for(i=0; i<h->nCHout; i++)
            free(h->buffer_syn[i]);
        free(h->buffer_win);
        free(h->win_sum);
        free(h->win_sum_cmplx_dummy);
        free(h->qmfTF_frame);
        free(h->qmfTF_frame_tmp);
        free(h->tmp_real_frame);
        free(h->tmp_imag_frame);

        /* For hybrid filtering */
        if(h->hybridmode){
            free(h->qmfDelayBuffer);
            free(h->hybBuffer);
            free(h->hybQmfTF_frame);
        }

        free(h);
        h=NULL;
        *phQMF = NULL;
    }
}

void qmf_analysis
(
    void * const hQMF,
    float** dataTD,
    int framesize,
    float_complex*** dataFD 
)
{
    qmf_data *h = (qmf_data*)(hQMF);
    int i, ch, t, nHops, band;
    float_complex subBands8[8], subBands2[2];
    const float_complex calpha = cmplxf(1.0f, 0.0f), cbeta = cmplxf(0.0f, 0.0f);

    saf_assert(framesize % h->hopsize == 0, "framesize must be multiple of hopsize");  
    nHops = framesize/h->hopsize;

    for(ch=0; ch<h->nCHin; ch++){
        for(t=0; t<nHops; t++){
            /* Shift samples to the right by one hopsize, and copy the current frame
               to the beginning */
            memmove(&(h->buffer_ana[ch][h->hopsize]), h->buffer_ana[ch], h->hopsize * 9 * sizeof(float));
            cblas_scopy(h->hopsize, &dataTD[ch][t*(h->hopsize)], -1, h->buffer_ana[ch], 1);

            /* Apply prototype filter/window */
            utility_svvmul(h->buffer_ana[ch], h->h_p, h->hopsize*10, h->buffer_win);

            /* Sum all 5 consecutive 1:2*hopsize */
            utility_svvadd(h->buffer_win, &(h->buffer_win[h->hopsize*2]), h->hopsize*2, h->win_sum);
            cblas_saxpy(h->hopsize*2, 1.0f, h->buffer_win + h->hopsize*4, 1, h->win_sum, 1);
            cblas_saxpy(h->hopsize*2, 1.0f, h->buffer_win + h->hopsize*6, 1, h->win_sum, 1);
            cblas_saxpy(h->hopsize*2, 1.0f, h->buffer_win + h->hopsize*8, 1, h->win_sum, 1);
            cblas_scopy(h->hopsize*2, h->win_sum, 1, h->win_sum_cmplx_dummy, 2);

            /* Apply complex-QMF analysis modulators  */
            cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, h->hopsize, 1, h->hopsize*2, &calpha,
                        FLATTEN2D(h->h_a), h->hopsize*2,
                        (float_complex*)h->win_sum_cmplx_dummy, 1, &cbeta,
                        h->qmfTF_frame, 1);

            /* Subdivide the lowest 3 bands */
            if(h->hybridmode){
                /* Shift buffer down by 1 frame */
                memmove(h->hybBuffer[ch][0], &(h->hybBuffer[ch][0][1]), (QMF_HYBRID_FILTER_LENGTH-1)*sizeof(float_complex));
                memmove(h->hybBuffer[ch][1], &(h->hybBuffer[ch][1][1]), (QMF_HYBRID_FILTER_LENGTH-1)*sizeof(float_complex));
                memmove(h->hybBuffer[ch][2], &(h->hybBuffer[ch][2][1]), (QMF_HYBRID_FILTER_LENGTH-1)*sizeof(float_complex));

                /* Append new frame to hybrid filtering buffer */
                h->hybBuffer[ch][0][QMF_HYBRID_FILTER_LENGTH-1] =  h->qmfTF_frame[0];
                h->hybBuffer[ch][1][QMF_HYBRID_FILTER_LENGTH-1] =  h->qmfTF_frame[1];
                h->hybBuffer[ch][2][QMF_HYBRID_FILTER_LENGTH-1] =  h->qmfTF_frame[2];

                /* Delay all the other QMF bands (i.e., the ones not being subdivided)
                 * so that they align with the hybrid bands in time: */
                for(i=0; i<h->hopsize - QMF_NBANDS_2_SUBDIVIDE; i++){
                    memmove(h->qmfDelayBuffer[ch][i], &(h->qmfDelayBuffer[ch][i][1]), ((QMF_HYBRID_FILTER_LENGTH-1)/2)*sizeof(float_complex));
                    h->qmfDelayBuffer[ch][i][(QMF_HYBRID_FILTER_LENGTH-1)/2] = h->qmfTF_frame[i+QMF_NBANDS_2_SUBDIVIDE];
                }

                /* Subdivide first QMF band into 8 subbands, and form hybrid bands 1-6 */
                cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, 8, 1, QMF_HYBRID_FILTER_LENGTH, &calpha,
                            h->fb8bandCoeffs, QMF_HYBRID_FILTER_LENGTH,
                            h->hybBuffer[ch][0], 1, &cbeta,
                            subBands8, 1);
                h->hybQmfTF_frame[0] = subBands8[6];
                h->hybQmfTF_frame[1] = subBands8[7];
                h->hybQmfTF_frame[2] = subBands8[0];
                h->hybQmfTF_frame[3] = subBands8[1];
#if _MSC_VER >= 1900
                h->hybQmfTF_frame[4] = ccaddf(subBands8[2], subBands8[5]);
                h->hybQmfTF_frame[5] = ccaddf(subBands8[3], subBands8[4]);
#else
                h->hybQmfTF_frame[4] = subBands8[2] + subBands8[5];
                h->hybQmfTF_frame[5] = subBands8[3] + subBands8[4];
#endif

                /* Subdivide second QMF band to get hybrid bands 7 and 8 */
                cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, 2, 1, QMF_HYBRID_FILTER_LENGTH, &calpha,
                            h->fb4bandCoeffs, QMF_HYBRID_FILTER_LENGTH,
                            h->hybBuffer[ch][1], 1, &cbeta,
                            subBands2, 1);
                h->hybQmfTF_frame[6] = subBands2[1]; /* Flipped! */
                h->hybQmfTF_frame[7] = subBands2[0];

                /* Subdivide third QMF band to get hybrid bands 9 and 10 */
                cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, 2, 1, QMF_HYBRID_FILTER_LENGTH, &calpha,
                            h->fb4bandCoeffs, QMF_HYBRID_FILTER_LENGTH,
                            h->hybBuffer[ch][2], 1, &cbeta,
                            subBands2, 1);
                h->hybQmfTF_frame[8] = subBands2[0];
                h->hybQmfTF_frame[9] = subBands2[1];

                /* The remaining bands are then just the delayed qmf bands, 4:end */
                cblas_ccopy(h->hopsize - QMF_NBANDS_2_SUBDIVIDE, FLATTEN2D(h->qmfDelayBuffer[ch]),
                            (QMF_HYBRID_FILTER_LENGTH-1)/2 + 1, &(h->hybQmfTF_frame[10]), 1);
            }

            /* copy to output */
            if(h->hybridmode){
                switch(h->format){
                    case QMF_BANDS_CH_TIME:
                        for(band=0; band<h->nBands; band++)
                            dataFD[band][ch][t] = h->hybQmfTF_frame[band];
                        break;
                    case QMF_TIME_CH_BANDS:
                        memcpy(dataFD[t][ch], h->hybQmfTF_frame, h->nBands*sizeof(float_complex));
                        break;
                }
            }
            else{
                switch(h->format){
                    case QMF_BANDS_CH_TIME:
                        for(band=0; band<h->nBands; band++)
                            dataFD[band][ch][t] = h->qmfTF_frame[band];
                        break;
                    case QMF_TIME_CH_BANDS:
                        memcpy(dataFD[t][ch], h->qmfTF_frame, h->nBands*sizeof(float_complex));
                        break;
                }
            }
        }
    }
}

void qmf_synthesis
(
    void * const hQMF,
    float_complex*** dataFD,
    int framesize,
    float** dataTD
)
{
    qmf_data *h = (qmf_data*)(hQMF);
    int ch, t, nHops, band;

    saf_assert(framesize % h->hopsize == 0, "framesize must be multiple of hopsize");
    nHops = framesize/h->hopsize;

    for(ch=0; ch<h->nCHout; ch++){
        for(t=0; t<nHops; t++){
            /* Load frequency domain data */
            if(h->hybridmode){
                switch(h->format){
                    case QMF_BANDS_CH_TIME:
                        for(band=0; band<h->nBands; band++)
                            h->hybQmfTF_frame[band] = dataFD[band][ch][t];
                        break;
                    case QMF_TIME_CH_BANDS:
                        memcpy(h->hybQmfTF_frame, dataFD[t][ch], h->nBands*sizeof(float_complex));
                        break;
                }

                /* Recombine the hybrid bands: */
#if _MSC_VER >= 1900
                h->qmfTF_frame[0] = ccaddf( h->hybQmfTF_frame[0], h->hybQmfTF_frame[1]);
                h->qmfTF_frame[0] = ccaddf( h->qmfTF_frame[0],    h->hybQmfTF_frame[2]);
                h->qmfTF_frame[0] = ccaddf( h->qmfTF_frame[0],    h->hybQmfTF_frame[3]);
                h->qmfTF_frame[0] = ccaddf( h->qmfTF_frame[0],    h->hybQmfTF_frame[4]);
                h->qmfTF_frame[0] = ccaddf( h->qmfTF_frame[0],    h->hybQmfTF_frame[5]);
                h->qmfTF_frame[1] = ccaddf( h->hybQmfTF_frame[6], h->hybQmfTF_frame[7]);
                h->qmfTF_frame[2] = ccaddf( h->hybQmfTF_frame[8], h->hybQmfTF_frame[9]);
#else
                h->qmfTF_frame[0] = h->hybQmfTF_frame[0]+h->hybQmfTF_frame[1]+h->hybQmfTF_frame[2]+
                                    h->hybQmfTF_frame[3]+h->hybQmfTF_frame[4]+h->hybQmfTF_frame[5];
                h->qmfTF_frame[1] = h->hybQmfTF_frame[6]+h->hybQmfTF_frame[7];
                h->qmfTF_frame[2] = h->hybQmfTF_frame[8]+h->hybQmfTF_frame[9];
#endif
                memmove(&(h->qmfTF_frame[3]), &(h->hybQmfTF_frame[10]), (h->hopsize - QMF_NBANDS_2_SUBDIVIDE)*sizeof(float_complex));
            }
            else{
                switch(h->format){
                    case QMF_BANDS_CH_TIME:
                        for(band=0; band<h->nBands; band++)
                            h->qmfTF_frame[band] = dataFD[band][ch][t];
                        break;
                    case QMF_TIME_CH_BANDS:
                        memcpy(h->qmfTF_frame, dataFD[t][ch], h->nBands*sizeof(float_complex));
                        break;
                }
            }

            /* Shift samples to the right by 2*hopsize */
            memmove(h->buffer_syn[ch] + h->hopsize*2, h->buffer_syn[ch], h->hopsize * 18 * sizeof(float));

            /* Apply complex-QMF synthesis modulators */
            cblas_scopy(h->hopsize, (float*)h->qmfTF_frame, 2, h->qmfTF_frame_tmp, 1); /* creal(h->qmfTF_frame) */
            cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, h->hopsize*2, 1, h->hopsize, 1.0f,
                        FLATTEN2D(h->h_s_real), h->hopsize,
                        h->qmfTF_frame_tmp, 1, 0.0f,
                        h->tmp_real_frame, 1);
            cblas_scopy(h->hopsize, &((float*)h->qmfTF_frame)[1], 2, h->qmfTF_frame_tmp, 1); /* cimag(h->qmfTF_frame) */
            cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, h->hopsize*2, 1, h->hopsize, 1.0f,
                        FLATTEN2D(h->h_s_imag), h->hopsize,
                        h->qmfTF_frame_tmp, 1, 0.0f,
                        h->tmp_imag_frame, 1);

            /* Append new synthesis frame */
            utility_svvsub(h->tmp_real_frame, h->tmp_imag_frame, h->hopsize*2, h->buffer_syn[ch]);

            /* Apply prototype filter/window */
            utility_svvmul(h->buffer_syn[ch], h->h_p, h->hopsize, h->buffer_win);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*3,  h->h_p + h->hopsize,   h->hopsize, h->buffer_win + h->hopsize);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*4,  h->h_p + h->hopsize*2, h->hopsize, h->buffer_win + h->hopsize*2);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*7,  h->h_p + h->hopsize*3, h->hopsize, h->buffer_win + h->hopsize*3);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*8,  h->h_p + h->hopsize*4, h->hopsize, h->buffer_win + h->hopsize*4);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*11, h->h_p + h->hopsize*5, h->hopsize, h->buffer_win + h->hopsize*5);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*12, h->h_p + h->hopsize*6, h->hopsize, h->buffer_win + h->hopsize*6);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*15, h->h_p + h->hopsize*7, h->hopsize, h->buffer_win + h->hopsize*7);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*16, h->h_p + h->hopsize*8, h->hopsize, h->buffer_win + h->hopsize*8);
            utility_svvmul(h->buffer_syn[ch] + h->hopsize*19, h->h_p + h->hopsize*9, h->hopsize, h->buffer_win + h->hopsize*9);

            /* Sum all 1:hopsizes to get output frame */
            utility_svvadd(h->buffer_win, h->buffer_win + h->hopsize,   h->hopsize, dataTD[ch] + t*(h->hopsize));
            cblas_saxpy(h->hopsize, 1.0f, h->buffer_win + h->hopsize*2, 1, dataTD[ch] + t*(h->hopsize), 1);
            cblas_saxpy(h->hopsize, 1.0f, h->buffer_win + h->hopsize*3, 1, dataTD[ch] + t*(h->hopsize), 1);
            cblas_saxpy(h->hopsize, 1.0f, h->buffer_win + h->hopsize*4, 1, dataTD[ch] + t*(h->hopsize), 1);
            cblas_saxpy(h->hopsize, 1.0f, h->buffer_win + h->hopsize*5, 1, dataTD[ch] + t*(h->hopsize), 1);
            cblas_saxpy(h->hopsize, 1.0f, h->buffer_win + h->hopsize*6, 1, dataTD[ch] + t*(h->hopsize), 1);
            cblas_saxpy(h->hopsize, 1.0f, h->buffer_win + h->hopsize*7, 1, dataTD[ch] + t*(h->hopsize), 1);
            cblas_saxpy(h->hopsize, 1.0f, h->buffer_win + h->hopsize*8, 1, dataTD[ch] + t*(h->hopsize), 1);
            cblas_saxpy(h->hopsize, 1.0f, h->buffer_win + h->hopsize*9, 1, dataTD[ch] + t*(h->hopsize), 1);
        }
    }
}

void qmf_channelChange
(
    void * const hQMF,
    int new_nCHin,
    int new_nCHout
)
{
    qmf_data *h = (qmf_data*)(hQMF);
    int i;

    if(h->nCHin!=new_nCHin){
        /* resize hybrid analysis buffers */
        if(h->hybridmode){
            h->qmfDelayBuffer = (float_complex***)realloc3d_r((void***)h->qmfDelayBuffer, new_nCHin, h->hopsize-QMF_NBANDS_2_SUBDIVIDE,
                                                              (QMF_HYBRID_FILTER_LENGTH-1)/2 + 1, h->nCHin, h->hopsize-QMF_NBANDS_2_SUBDIVIDE,
                                                              (QMF_HYBRID_FILTER_LENGTH-1)/2 + 1, sizeof(float_complex));
            h->hybBuffer = (float_complex***)realloc3d_r((void***)h->hybBuffer, new_nCHin, QMF_NBANDS_2_SUBDIVIDE,
                                                         QMF_HYBRID_FILTER_LENGTH, h->nCHin, QMF_NBANDS_2_SUBDIVIDE,
                                                         QMF_HYBRID_FILTER_LENGTH, sizeof(float_complex));

            /* zero any new channels */
            for(i=h->nCHin; i<new_nCHin; i++){
                memset(FLATTEN2D(h->qmfDelayBuffer[i]), 0, (h->hopsize-QMF_NBANDS_2_SUBDIVIDE) * ((QMF_HYBRID_FILTER_LENGTH-1)/2 + 1) * sizeof(float_complex));
                memset(FLATTEN2D(h->hybBuffer[i]), 0, QMF_NBANDS_2_SUBDIVIDE * QMF_HYBRID_FILTER_LENGTH * sizeof(float_complex));
            }
        }

        /* resize analysis buffers */
        for(i=new_nCHin; i<h->nCHin; i++)
            free(h->buffer_ana[i]);
        h->buffer_ana = (float**)realloc1d(h->buffer_ana, sizeof(float*)*new_nCHin);
        for(i=h->nCHin; i<new_nCHin; i++)
            h->buffer_ana[i] = (float*)calloc1d(h->hopsize * 10,sizeof(float));

        h->nCHin = new_nCHin;
    }

    if(h->nCHout!=new_nCHout){
        /* resize synthesis buffers */
        for(i=new_nCHout; i<h->nCHout; i++)
            free(h->buffer_syn[i]);
        h->buffer_syn = (float**)realloc1d(h->buffer_syn, sizeof(float*)*new_nCHout);
        for(i=h->nCHout; i<new_nCHout; i++)
            h->buffer_syn[i] = (float*)calloc1d(h->hopsize * 20, sizeof(float));

        h->nCHout = new_nCHout;
    }
}

void qmf_clearBuffers
(
    void * const hQMF
)
{
    qmf_data *h = (qmf_data*)(hQMF);
    int i;

    /* flush analysis buffers */
    for(i=0; i<h->nCHin; i++){
        memset(h->buffer_ana[i], 0, h->hopsize * 10 * sizeof(float));
        if(h->hybridmode){
            memset(FLATTEN3D(h->qmfDelayBuffer), 0, h->nCHin*(h->hopsize-QMF_NBANDS_2_SUBDIVIDE)* ((QMF_HYBRID_FILTER_LENGTH-1)/2 + 1)*sizeof(float_complex));
            memset(FLATTEN3D(h->hybBuffer), 0, h->nCHin*QMF_NBANDS_2_SUBDIVIDE*QMF_HYBRID_FILTER_LENGTH*sizeof(float_complex));
        }
    }

    /* flush synthesis buffers */
    for(i=0; i<h->nCHout; i++)
        memset(h->buffer_syn[i], 0, h->hopsize * 20 * sizeof(float));
}

int qmf_getProcDelay
(
    void * const hQMF
)
{
    qmf_data *h = (qmf_data*)(hQMF);
    return h->procDelay;
}

int qmf_getNBands
(
    void * const hQMF
)
{
    qmf_data *h = (qmf_data*)(hQMF);
    return h->nBands;
}

void qmf_getCentreFreqs
(
    void * const hQMF,
    float fs,
    int nBands,
    float* centreFreq
)
{
    qmf_data *h = (qmf_data*)(hQMF);
    int i, j, hopsize;
    float* cutoffs;
    float centreFreqs_qmf[QMF_NBANDS_2_SUBDIVIDE];

    saf_assert(nBands==h->nBands, "Just to check that \"centreFreq\" is of correct length");  

    /* QMF cutoff frequencies */
    hopsize = h->hopsize;
    cutoffs = malloc1d((hopsize+1)*sizeof(float));
    for(i=0; i<hopsize+1; i++)
        cutoffs[i] = i* fs/(2.0f*(float)hopsize);

    if(h->hybridmode){
        /* Convert the first 3 qmf centre frequencies to the 10 hybrid centre frequencies */
        for(i=0; i<QMF_NBANDS_2_SUBDIVIDE; i++)
            centreFreqs_qmf[i] = cutoffs[i+1] - (cutoffs[i+1]-cutoffs[i])/2.0f;
        cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, 10, 1, QMF_NBANDS_2_SUBDIVIDE, 1.0f,
                    (const float*)__qmf2hybCentreFreq, QMF_NBANDS_2_SUBDIVIDE,
                    centreFreqs_qmf, 1, 0.0f,
                    centreFreq, 1);

        /* Remaining centre frequencies are then QMF centre frequencies 4:end  */
        for(i=10, j=QMF_NBANDS_2_SUBDIVIDE; i<nBands; i++, j++)
            centreFreq[i] = cutoffs[j+1] - (cutoffs[j+1]-cutoffs[j])/2.0f;
    }
    else{
        /* Centre frequencies are defined as directly inbetween the cutoff frequencies */
        for(i=0; i<hopsize; i++)
            centreFreq[i] = cutoffs[i+1] - (cutoffs[i+1]-cutoffs[i])/2.0f;
    }

    free(cutoffs);
}

void qmf_FIRtoFilterbankCoeffs
(
    float* hIR /*N_dirs x nCH x ir_len*/,
    int N_dirs,
    int nCH,
    int ir_len,
    int hopSize,
    int hybridmode,
    float_complex* hFB /* nBands x nCH x N_dirs */
)
{
    int i, j, t, nd, nm, nTimeSlots, ir_pad, nBands;
    int* maxIdx;
    float maxVal, idxDel, irFB_energy, irFB_gain, phase;
    float* centerImpulse, *centerImpulseFB_energy, *ir;
    float_complex cross;
    float_complex* centerImpulseFB, *irFB;

    nBands = hopSize + (hybridmode ? 7 : 0);
    ir_pad = 1024;//+512;
    nTimeSlots = (SAF_MAX(ir_len,hopSize)+ir_pad)/hopSize;
    maxIdx = calloc1d(nCH,sizeof(int));
    centerImpulse = calloc1d(SAF_MAX(ir_len,hopSize)+ir_pad, sizeof(float));

    /* pick a direction to estimate the center of FIR delays */
    for(j=0; j<nCH; j++){
        maxVal = 2.23e-13f;
        for(i=0; i<ir_len; i++){
            if(hIR[j*ir_len + i] > maxVal){
                maxVal = hIR[j*ir_len + i];
                maxIdx[j] = i;
            }
        }
    }
    idxDel = 0.0f;
    for(j=0; j<nCH; j++)
        idxDel += (float)maxIdx[j];
    idxDel /= (float)nCH;
    idxDel = (idxDel + 1.5f);

    /* ideal impulse at mean delay */
    centerImpulse[(int)idxDel] = 1.0f;

    /* analyse impulse with the filterbank */
    centerImpulseFB = malloc1d(nBands*nTimeSlots*nCH*sizeof(float_complex));
    qmfAnalyse(centerImpulse, SAF_MAX(ir_len,hopSize)+ir_pad, 1, hopSize, hybridmode, centerImpulseFB);
    centerImpulseFB_energy = calloc1d(nBands, sizeof(float));
    for(i=0; i<nBands; i++)
        for(t=0; t<nTimeSlots; t++)
            centerImpulseFB_energy[i] += powf(cabsf(centerImpulseFB[i*nTimeSlots + t]), 2.0f);

    /* initialise FB coefficients */
    ir = calloc1d( (SAF_MAX(ir_len,hopSize)+ir_pad) * nCH, sizeof(float));
    irFB = malloc1d(nBands*nCH*nTimeSlots*sizeof(float_complex));
    for(nd=0; nd<N_dirs; nd++){
        for(j=0; j<ir_len; j++)
            for(i=0; i<nCH; i++)
                ir[j*nCH+i] = hIR[nd*nCH*ir_len + i*ir_len + j];
        qmfAnalyse(ir, SAF_MAX(ir_len,hopSize)+ir_pad, nCH, hopSize, hybridmode, irFB);
        for(nm=0; nm<nCH; nm++){
            for(i=0; i<nBands; i++){
                irFB_energy = 0;
                for(t=0; t<nTimeSlots; t++)
                    irFB_energy += powf(cabsf(irFB[i*nTimeSlots*nCH + t*nCH + nm]), 2.0f); /* out_nBands x nTimeslots x nCH */
                irFB_gain = sqrtf(irFB_energy/SAF_MAX(centerImpulseFB_energy[i], 2.23e-8f));
                cross = cmplxf(0.0f,0.0f);
                for(t=0; t<nTimeSlots; t++)
                    cross = ccaddf(cross, ccmulf(irFB[i*nTimeSlots*nCH + t*nCH + nm], conjf(centerImpulseFB[i*nTimeSlots + t])));
                phase = atan2f(cimagf(cross), crealf(cross));
                hFB[i*nCH*N_dirs + nm*N_dirs + nd] = crmulf( cexpf(cmplxf(0.0f, phase)), irFB_gain);
            }
        }
    }

    /* clean-up */
    free(maxIdx);
    free(centerImpulse);
    free(centerImpulseFB_energy);
    free(centerImpulseFB);
    free(ir);
    free(irFB);
}
